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It was recently proved that solitons embedded in the spectrum of linear waves may exist in discrete 
systems, and explicit solutions for isolated unstable embedded lattice solitons (ELS) of a differential- 
difference version of a higher-order NLS equation were found [Physica D 197 (2004) 86]. The 
discovery of these ELS gives rise to relevant questions such as the following: are there continuous 
families of ELS?, can ELS be stable?, is it possible for ELS to move along the lattice?, how do ELS 
interact?. The present work addresses these questions by showing that a novel equation (a discrete 
version of a complex modified KdV equation which includes next-nearest-neighbor couplings) has 
a two-parameter continuous family of exact ELS. These solitons can move with arbitrary velocities 
across the lattice, and the numerical simulations demonstrate that these ELS are completely stable. 
Moreover, the numerical tests show that these ELS are robust enough to withstand collisions, and 
the result of a collision is only a shift in the positions of the solitons. The model may apply to 
the description of a Bose-Einstein condensate with dipole-dipole interactions between the atoms, 
trapped in a deep optical-lattice potential. 
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PACS numbers: 05.45.Yv, 63.20.Ry, 45.05. +x, 05.45.-a 

In nonlinear systems where solitons can exist, 
the propagation of small-amplitude linear waves, 
which obey the linearized version of the nonlin- 
ear equations, is possible too. However, for a 
soliton to exist, it is absolutely necessary that no 
resonances occur between the soliton and these 
linear waves. Otherwise, the soliton would de- 
cay due to an energy transfer towards the linear 
waves. Based on this no-resonance argument, it 
was frequently assumed that the solitons' internal 
frequencies could not be contained in the linear 
spectrum of the system, i.e., they could not lie 
within the band of frequencies permitted to linear 
waves. However, at the end of the nineties excep- 
tions to this rule were found, and a special type 
of solitons were discovered, which do not resonate 
with linear waves, in spite of having frequencies 
immersed in the spectrum of these waves. In 1999 
these peculiar solitary waves were given the name 
of embedded solitons (ES), and in the following 
years a number of models supporting ES were 
identified. Most of these models describe continu- 
ous systems. However, some examples of discrete 
ES were recently found too. These embedded lat- 
tice solitons (ELS) are isolated solutions which 
are stable against small perturbations in the lin- 
ear approximation, but are nonlinearly unstable. 
The discovery of these isolated unstable ELS trig- 
gered the search for models admitting continuous 
families of stable ELS. In this article we present 
a novel differential-difference equation which has 
a two-parameter continuous family of exact ELS. 
These solitons can move with arbitrary velocities 
across the lattice, and the numerical simulations 



show that they are completely stable solutions. 



I. INTRODUCTION 

In the late nineties, a new type of solitary waves was 
identified and given the name of embedded solitons (ES). 
At first, ES were found in nonlinear-optical models in- 
cluding quadratic (x'^-*) nonlinearities [iMqIj and later 
they were discovered in hydrodynamic [3 and liquid- 
crystal dynamical models [g- Quite recently, ES were 
found in discrete systems too, viz., in a finite-difference 
version of a higher-order NLS (nonlinear Schrodinger) 
equation and in a model of an array of linearly- 
coupled w aveg uides with the x''^'' and x^'^'' (cubic) non- 
linearities ■ The distinctive feature of the ES is that 
their internal frequencies and, in some cases, also their 
velocities, fall in a spectral region occupied by linear wave 
modes, where regular solitons cannot exist, as they would 
decay into radiation. That is, ES are exceptional states 
that do not emit linear waves, despite being in resonance 
with them. A qualitative explanation to the existence of 
ES was proposed in Ref. j3| (see also 0), and more 
quantitative arguments were provided in Ref. [T^ . 

A number of generic features of the ES found in contin- 
uous systems have been identified. In particular, ES are 
usually isolated solutions, although continuous families 
of ES may exist too, under special symmetry conditions 
0,0, 13- The ES may be stable against small perturba- 
tions in the linear approximation, but general arguments 
and direct simulations demonstrate that they are usually 
unstable at the second order, which is why they are often 
called "semi-stable" solutions. Nevertheless, examples of 
completely stable ES have been also found [MEiEli and 
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their stability has been demonstrated in Refs. [l^lTsll. 

Much less is known about embedded lattice solitons 
(ELS) in discrete systems. Only a few examples of models 
admitting isolated nonlinearly unstable ELS have been 
found 9j .1Q | . Some issues which remain to be analyzed 
are, for example, the possibility of the existence of con- 
tinuous families of ELS, and to determine if they may be 
completely stable. Furthermore, a challenging question 
is whether ELS may move across the lattice with arbi- 
trary constant velocities. We address these issues in the 
present work, introducing a new nonlinear lattice model 
which gives rise to exact solutions for ELS which feature 
all these properties. 

The paper is organized as follows. In Sec. 2 we present 
the new model and its continuous family of moving lattice 
solitons. In Sec. 3 we show that these moving solitons 
are indeed ELS. In Sec. 4 the behavior of perturbed and 
colliding ELS is investigated by means of direct simula- 
tions which, in particular, demonstrate that more gen- 
eral solitons with different amplitudes exist too. Section 
5 concludes the paper. 



II. THE MODEL AND ITS FAMILY OF 
MOVING LATTICE-SOLITON SOLUTIONS 

Motivated by the recent discoveries of continuous fam- 
ilies of ES in the complex modified Korteweg-de Vries 
(cmKdV) equation '3| and in a generalized third-order 
nonlinear Schrodinger equation |l3l ITsI l , in the present 
paper we will consider a differential-difference equation 
of the following form: 

f„ -2sq {rn+i - rn-i) + eo {rn+2 ~ rn-2) 

(1) 

where r„ (t) is a complex variable, the overdot stands 
for d/dt, and Eq, p and S are real constants. At fixed t, 
the discrete values r„ (t) can be considered as the val- 
ues of a continuous function u (x, t), evaluated at equally 
spaced lattice points. In other words, we can consider 
that r„ (t) = u{x — uAx, t), where Aa; is the lattice spac- 
ing. If we now consider the limit when Ax ^ 0, a stan- 
dard limiting procedure (described in detail in Ref. 0) 
shows that the Eq. goes over into the equation: 

where e = 2eo i^x)^ and 7 = 12Aj:. If p — 25 = 6 
Eq. (0) is precisely the cmKdV equation studied in 0. 
It is worth observing that the coefficient in front of the 
nonlinear term in the cmKdV equation may be positive 
or negative, depending on the physical application. For 
example, in the description of light propagation in liquid- 
crystal waveguides, this coefficient is positive 8]. On the 
other hand, the derivation of the cmKdV equation for 



strongly dispersive waves in a weakly nonlinear medium 
by means of a multiple-scale expansion, presented in Ref. 
[l^ . leads to a negative nonlinear coefficient. 

If we neglect the nonlinear terms in Eq. we arrive 
at a linear equation which describes a dispersive lattice 
{i.e., a discrete version of the linear dispersive equation 
Ut = EUxxx)- However, for the existence of ELS solutions, 
both nonlinear terms in Eq. are needed. 

By means of the staggering transformation, 

r„(i) = (-*)"i?„(t) (3) 
Eq. can be cast in an alternative form, 

iRn + 2eo {Rn+l + Rn-l) ^ P \Rn'^ {Rji+1 + Rn-l) 
+ ieo{Rn+2- Rn-2) -l5\Rn\''{Rn+2- Rn-2) -=^0. (4) 

In this form, it is an extension of the well-known inte- 
grable Ablowitz-Ladik (AL) model that corresponds 
to the terms in the first row in Eq. Q). The first (lin- 
ear) additional term accounts for the higher- order lattice 
diffraction 0|. As shown in Ref. 0, the third-order 
lattice diffraction can be induced by biasing a lattice 
(with the coupling between next-nearest neighbors) by a 
phase gradient applied along it (in particular, by means 
of the diffraction-management technique in an array of 
nonlinear optical waveguides Tflj). Then, the nonlinear 
terms in Eq. Q may be regarded as generic nonlinear 
corrections to the second- and third-order lattice diffrac- 
tion. In physical systems, nonlinear terms of these types 
should be induced by long-range nonlinear interactions 
between sites of the lattice. A realization of this is known 
in Bose-Einstein condensates (BECs) with dipole-dipole 
interactions between atoms, trapped in a deep periodic 
potential (which is created by means of an optical lattice) 
|20l | . Note that Eq. does not include any on-site non- 
linearity, which can be eliminated in the BEG by means 
of the known technique of the Feshbach resonance 21] . 

A basic finding of this work is that Eq. Q has a 
family of exact soliton solutions, provided that e^p < 
and S > 2p, of the following form: 

r„ (t) = AsolSech (^^^^^^) ^^«It + nQA.)^ (5) 

where Q and Aa; are free parameters, and the values of 
Asoi, w, a and q are determined by the following equa- 
tions: 

AL = '-^(2--), (6) 
P \ PJ 

— = -2£oCos(2QAa;)— smh(^— j 

-F4eocos (QAx) -j— sinh , (8) 

Aa; \ w J 
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q = 2eosin (2QAx) cosh ) 

\ w J 

f Ax 

— 4£osin (QAx) cosh 

\ w 



(9) 



For fixed vahies of the coefficients Sq, p and S, the Eqs. 
©-(PI define a continuous family of moving lattice soli- 
tons, all with the same amplitude, due to Eq. ©, but 
different widths w and velocities V ~ a/ Ax (the lat- 
ter being defined with respect to the discrete coordinate 
n). Because the value of Ax/w is fixed by Eq. ||7J), and 
|cos {2QAx)\ < 1, Eq. (jHJ shows that the velocity always 
takes values within a limited interval, 



a 

Ax 



„, , w . , f2Ax 

<2|eo|-^smh 

Ax \ w 



(10) 



It is relevant to stress that, contrary to a freque ntly 
made assumption, explicit counter-examples |24l l2q| 
demonstrate that the existence of exact solutions for 
moving solitons does not imply the integrability of the 
corresponding nonlinear lattice system. Nevertheless, the 
availability of exact soliton solutions capable of moving 
with arbitrary velocities across the lattice is a highly non- 
trivial feature, which permits soliton collisions to occur. 
This is an interesting process which shall be studied in 
Sec. 4. 

To close this section we would like to notice that in 
the particular case when p = 26 the Eq. ^ can be de- 
rived from a Hamiltonian function (which is a dynamical 
invariant of the equation), provided that suitable non- 
canonical Poisson brackets are introduced. As in this 
case {p = 26) no soliton solutions for Eq. ^ have been 
found, we will not analyze this case here. However, the 
existence of a Hamiltonian structure for Eq. ^ is an 
interesting issue which merit further studies. 



III. SOLITON EMBEDDING 

We have seen that the Eq. has a continuous fam- 
ily of moving lattice solitons. Now we aim to determine 
whether these solitons are embedded in the spectrum of 
the radiation modes. To this end it is necessary to check 
if the intrinsic frequency of the soliton ((S} is contained 
(embedded) in the band of frequencies permitted to the 
small-amplitude linear waves (phonon modes) which are 
able to propagate in the lattice. If the soliton's frequency 
falls into this band, it will be important to determine 
the wavenumbers corresponding to all the phonon modes 
which have the same frequency as the soliton, since these 
modes will be resonantly excited if the soliton is per- 
turbed. 

In order to determine the correct resonant wavenum- 
bers it is necessary to take into account that the soli- 
tons defined by Eqs. ©-lini) are moving solitons. Due to 
this movement, the soliton's intrinsic frequency must be 
compared to the linear spectrum of the phonon waves in 



the reference frame which moves along with the soliton 
(the transition to the moving reference frame, which is 
generally impossible in a nonlinear lattice model, is al- 
ways possible in the linearized equation for the phonon 
modes). In this moving reference frame the coordinate 
nAx is replaced by nAx — at, thus implying that the 
accordingly defined frequency includes the Doppler shift. 
Additional arguments (based on the Fourier expansion of 
the fields) showing that the resonance between moving 
solitons and linear waves must be studied in the moving 
reference frame were given in Refs. |l2j and B and, in a 
more mathematically rigorous form, in Ref. [26j. 

According to the previous paragraph, we 
look for phonon modes in the form: r„ (t) = 
^ exp[int — ik {nAx — at)], with ^ being an arbitrary 
small amplitude. The substitution of this expression in 
the linearized form of Eq. Q leads to the dispersion 
relation: 



n (k) = 4eo sin {kAx) [1 — cos (kAx)] — ka, 



(11) 



the last term being the above-mentioned Doppler shift. 
The moving soliton will be embedded if its intrinsic fre- 
quency, which is — (q -|- aQ) according to Eq. |SJl, is con- 
tained in the range of values of the function —Q (fc) with 
the real argument k varying within — oo < k < +oo. 
Since this range is the entire real axis (if a 0), all 
the moving lattice solitons given by Eq. Q are indeed 
embedded sohtons. 

In the particular case of a = 0, the corresponding qui- 
escent soliton will be embedded if there is a real solution 
of the corresponding resonance condition, which follows 
from equating expressions © and l(TT|) : 



q — 4eo sin (kAx) [I — cos (kAx)] 



(12) 



It is easy to see that real solutions to Eq. p2|) exist for 
l^l < 3\/3|£o|, otherwise the static soliton will not be 
embedded. 

We would like to observe that some of the embedded 
solitons of the family defined by the Eqs. (0-® move 
with velocities which are contained within the range of 
phase velocities permitted to the phonon modes. From 
(|11|1 it follows that the phase velocities of these modes 
(defined in the laboratory reference frame) are defined 
by: 



vik) 



4£o 

kAx 



sin (kAx) [1 - cos (kAx)] . (13) 



Consequently, the soliton's velocity {V = a/ Ax) will be 
within the range of this function if the equation: 



4£o 
k 



sin {kAx) [1 — cos {kAa 



(14) 



has real solution for k. Depending on the values of Eq, p, 
6, Q and Ax, real solutions to Eq. H14|) may or may not 
exist. 
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From the previous paragraph it follows that the soli- 
tons of Eq. ^ can be classified in two groups. The first 
one contains those solitons whose internal frequencies lie 
within the range of the function — fl (fc), and whose ve- 
locities are not contained in the range of the function 
The second group includes those solitons which, 
in addition of having wavenumbers within the range of 
—Q (fc), have velocities which are contained in the range 
of the function (|13|l . In the following, the solitons of the 
first group will be called single- embedded solitons, and 
the solitons of the second group will be referred to as 
double- embedded solitons, to emphasize that in this case 
both, the wavenumbers and the velocities of the solitons, 
are contained in the corresponding linear spectra. 



IV. PERTURBED AND COLLIDING LATTICE 
SOLITONS 

A. Spectrum of emission from perturbed 
embedded solitons 

Since the intrinsic frequencies of the soliton solutions 
of Eq. 1^ fall within the phonon spectrum of the sys- 
tem, one might expect that resonances should occur be- 
tween these solitons and the phonon modes, giving rise 
to the emission of resonant radiation. However, the dis- 
tinctive feature of the ES is precisely the absence of this 
radiation. In the case of Eq. ||2Jl, the mechanism that 
suppresses the resonance between the ES and the linear 
waves can be understood if one calculates the Fourier 
transform u{k,uj) of the field u{x,t). As demonstrated 
in Ref. the calculation shows that u (fc, w) is propor- 
tional to a quotient, P2 (w) /Pi (w), of two polynomials. 
Pi (w) is an image of the linear part of the equation, and 
resonances occur at frequencies for which Pi (w) vanishes, 
while P2 (w) depends on the nonlinear part of the equa- 
tion. In order to calculate the Fourier transform of the 
nonlinear term it is considered that u{x,t) is a single- 
humped function (such as a hyperbolic secant times a 
complex exponential). In this way an explicit form for 
P2 (w) can be obtained (see Eq. (52) in Ref. 0). The 
crucial point which explains why the ES do not radiate is 
that P2 (w) vanishes at the same points where Pi (w) = 0, 
provided that u{x,t) is exactly an ES, and consequently 
the singularities in the expression P2 (ui) /Pi (uj), which 
lead to the resonances, cancel out. However, if u{x,t) 
is an altered (perturbed) ES, the exact cancelation does 
not take place, and hence transient emission from the 
perturbed ES will take place at the frequencies for which 
Pi (oj) = 0. Generally, if the soliton's intrinsic frequency 
is uj = — {q -\- aQ) [as in Eq. Q], and (fc) is the linear 
dispersion relation for the phonon modes in the reference 
frame moving along with the soliton, such as in Eq. (|ll|l . 
then Pi (w) = q-\- aQ — Vlik). Consequently, a perturbed 
ES is expected to emit radiation at wavenumbers —k de- 
termined by the resonance condition (fc) = 5-1- aQ. In 
other words, the resonant wavenumbers are defined by 



the equation: 

q -\- aQ — 4eosin (kAx) [1 — cos {kAx)] — ka. (15) 

To confirm that the perturbed solitons emit radiation 
at these wavenumbers, let us analyze two examples. In 
the first case let us consider the following parameters: 
£0 = 1, p = -0.4, 5 = -0.85, Q = 0.61 and Ax = 3. 
In this case, the soliton ||SJ) is a double-embedded one, 
with Asoi = 0.56, w = 12.12, a = 7.75 and q = -5.1. 
To verify that this soliton is double-embedded, in Fig. 1 
we can see the graph of the function v (k) Ax [with v (k) 
given by Eq. H13|) ]. Since the maximum of this function 
is 7.967, the soliton's parameter a = 7.75 lies within the 
range of v (k) Ax, thus implying that the soliton is indeed 
double-embedded. To perturb this soliton, we solved Eq. 
^ numerically with the initial condition: 

r„ (t = 0) = Ao scch j e*"Q^^ (16) 

where the initial amplitude, Aq — 0.64, is larger than the 
amplitude of the exact soliton, A^oi = 0.56. According 
to Eq. (HSl, the perturbed soliton will emit radiation at 
the wavenumbers: 

= -0.114, -0.073, -0.008, +0.103, +0.129. 

The numerical solution of the Eq. Q with the initial 
condition H16|l is shown in Fig. 2 (for t — 20), and the 
corresponding power spectrum, |r(/c)|^, is displayed in 
Fig. 3 within the Brillouin zone, — tt/Ax < k < -\-tt/Ax, 
where the Fourier transform of the lattice field is defined 
as: 

r(fc) = r„ exp (iknAx) . (18) 

n 

The spectrum shown in Fig. 3 contains a central com- 
ponent, which is essentially the Fourier transform of 
the unperturbed soliton, and two symmetric radiation 
bands, approximately located at 0.07 < |fc/27r| < 0.13. 
These radiation bands contain all the expected resonant 
wavenumbers, except the value —0.008, which falls within 
the central band. The symmetry of the radiation bands is 
related to the fact that the oscillating behavior (in time) 
of the soliton solutions of Eq. ^ is described by a com- 
plex exponential exp (iujot), whose real part oscillates as 
cos(±woi)- Due to this fact, a perturbed soliton is able 
to excite the phonon modes with frequencies ±ujo, thus 
accounting for the occurrence of symmetrically refiected 
resonance numbers. 

As a second example, let us perturb the single- 
embedded soliton corresponding to the parameters: £0 — 
-4.8, p = 4.8, 6 = 10.3, Q = 7.9 and Ax = 0.4. In 
this case the parameters of the soliton are: ^soi = 0.381, 
w = 1.503, a = 15.812 and q = -0.771. To confirm that 
this soliton is single-embedded, in Fig. 4 we display the 
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phase velocity (times Ax) of the radiation waves as given 
by Eq. H13I) . From this figure it is evident that the soH- 
ton's parameter a = 15.812 Hes outside the range of the 
function v (k) Ax. 

The resonance condition (|15|l defines in this case only 
one resonant wavenumber located at —k/2'K = —1.249, 
which is very close to one of the edges of the Brillouin 
zone |fc/27r| < l/(2Ax) = 1.25. As in the previous ex- 
ample, we also expect a symmetrically reflected radiation 
peak located at —kjl-n = +1.249, due to a complex fac- 
tor describing the oscillating behavior of the perturbed 
soliton. 

The numerical solution of Eq. |^ corresponding to an 
initial condition of the form with Aq = 0.35 (which 
is smaller than the amplitude of the exact soliton) shows 
that the perturbed pulse emits a radiation wavetrain to 
the left, as seen in Fig. 5. The spectrum of the numer- 
ical solution, shown in Fig. 6, exhibits two symmetric 
radiation peaks at — fc/27r = ±1.23, which are very close 
to the expected values. 



B. Robustness of the embedded solitons 

The emission of radiation from a perturbed linearly 
stable ELS may result in its full decay, and consequently 
its stability (or, better said, robustness) must be tested 
in long simulations. It is relevant to mention that in the 
continuum limit, the ES solutions of Eq. Q are stable 
because this equation gives rise to a continuous family 
of ES with arbitrary energy, and hence a perturbed soli- 
ton can shed off a part of its energy and settle down to 
an ES solution with a smaller energy In continuum 
models supporting isolated ES, the perturbed pulse may 
not be able to find a stationary soliton with a lower en- 
ergy toward which it could relax, and hence it may go 
on radiating till complete decay. This scenario, referred 
to as the "semi-stability" (alias nonlinear instability) of 
the ES, really occurs in the continuum model with the 
combined x^"^^ ■ X^'^^ nonlinearity [3LflT|. 

In the present case, the stability of the ELS may be 
secured if a perturbed soliton can find a way to relax, 
after shedding off some radiation, to a new stationary 
soliton shape. Within the bounds of the exact solution 
family given by Eqs. 10-®, this process does not seem 
plausible, as all the exact solutions have a unique value 
of the amplitude. Nevertheless, it may happen that the 
exact family is only a subset of a more general manifold 
of soliton solutions. The possible existence of such an 
extended soliton family would not only provide for the 
stability as explained above, but it is a fundamentally 
important issue by itself. 

To clarify the eventual fate of the perturbed solitons 
of Eq. Ip^lj we now display some typical examples. In the 
first place, let us consider the double-embedded soliton 
corresponding to Eq — I, p — —0.4, 5 — —0.85, Q = 0.61 
and Ax = 3. As stated before, the shape of this soliton 
is defined by the parameters Agoi = 0.56, w = 12.12 and 



a = 7.75. If this soliton is perturbed by increasing its 
amplitude [for instance, taking the initial condition l|16|) 
with Aq = 0.64], the perturbed pulse starts to oscillate. 
In the course of the evolution, the oscillations gradually 
fade away due to the radiation loss, as shown by the 
upper curve in Fig. 7, and the perturbed soliton tends 
to settle down into a stationary shape, which is quite 
possible, as the initial perturbation increased its norm. 

On the other hand, if the height of the initial pulse 
is smaller than the amplitude of the exact soliton (for 
example, Aq — 0.48), the norm of the perturbed pulse 
becomes smaller than the exact soliton's norm, and there- 
fore the pulse will not be able to relax back to the exact 
soliton. As explained above, such a perturbation might 
lead to the destruction of the soliton, unless it can find 
another stationary shape to which it may relax to. The 
evolution of this perturbed soliton was numerically simu- 
lated, and the lower curve in Fig. 7 shows the behavior of 
the soliton's amplitude. As we can see in this figure, the 
perturbed soliton does not decay indefinitely. It again 
exhibits a damped oscillatory behavior (although the os- 
cillation period is much longer than in the previous case, 
when the perturbation increased the norm), and it finally 
approaches a new equilibrium state. The behavior exhib- 
ited in Fig. 7 indicates that the double-embedded lattice 
solitons are completely stable pulses and also, as conjec- 
tured above, that the exact solitons seem to be just a 
subset of a broader family of discrete solitons. The com- 
plete characterization of this broader family lies outside 
the scope of the present work and it will not be addressed 
here. 

The stability of the single-embedded lattice solitons of 
Eq. J^l is even clearer. Let us consider, for instance, the 
soliton corresponding to Sq — —4.8, p ~ 4.8, 6 = 10.3, 
Q = 7.9, and Ax = 0.4. In this case, A^oX = 0.381, 
w = 1.503, q = -0.771 and a = 15.812. As the value 
of the parameter a lies outside the range of the func- 
tion u(fc)Ax (see Fig. 4), the soliton is indeed single- 
embedded (i.e., only in terms of its frequency, but not ac- 
cording to its velocity). To perturb this soliton, we first 
took the initial condition (|16|l with an increased ampli- 
tude, Aq — 0.415. This perturbed pulse stabilizes itself 
quite rapidly, as we can see from the time evolution of 
its height, which is presented in the upper curve of Fig. 
8. In a similar way, if the initial amplitude is decreased 
to Aq = 0.345 (recall that the amplitude of the exact 
soliton is Aso\ = 0.381), the pulse also relaxes quickly, its 
height evolving as shown by the lower curve in Fig. 8. 
This figure indicates that the single-embedded soliton of 
Eq. tested in these simulations is a stable solution. 
Figure 8 also shows that when the initial pulse is higher 
(lower) than the exact soliton, the amplitude of the es- 
tablished steady-state pulse will be still larger (smaller) 
than the amplitude ^soi of the exact soliton solution. A 
similar behavior was observed in the case of the ES of 
Eq. and an explanation of this behavior (by means 
of a variational analysis) was provided in Ref. [Sj. 

A comparison of Figs. 7 and 8 reveals an interesting 
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difference: while the amphtude of the perturbed double- 
embedded soliton shown in Fig. 7 exhibits well-defined 
oscillations, the amplitude of the single-embedded soliton 
seen in Fig. 8 evolves almost monotonically. This quali- 
tative difference may be related to the fact that the reso- 
nance condition H15|l determines five resonant wavenum- 
bers in the case studied in Fig. 7, and only one in the 
case corresponding to Fig. 8. Actually, a precise expla- 
nation of the different behaviors exhibited in Figs. 7 and 
8 remains a challenge for future work. 

We again stress that the perturbed lattice solitons of 
Eq. stabilize as new stable pulses whose heights are 
different from the unique value of the amplitude of the 
exact soliton solutions given by Eq. ©. As said above, 
this observation strongly suggests that the exact solu- 
tions given by Eqs. ©-© constitute just a particular 
subset of a wider soliton family of the Eq. (QJ, which is 
a subject for a separate work. 

Our simulations, performed at many values of the pa- 
rameters, have never revealed an unstable soliton, even if 
the initial perturbation was actually large. An additional 
example is displayed in Fig. 9, for Sq = —1, p = 0.4, 
S = 0.85, Ax — 3 and Q = 0.61. The exact soliton cor- 
responding to these values has Agoi = 0.56, w = 12.12, 
a = 7.75, and q = —5.103, while the initial condition was 
taken in the form of Eq. H16|l with Aq = 0.48, which is 
essentially lower than the amplitude of the exact solu- 
tion. In this case again, the amplitude of this perturbed 
pulse exhibits damped oscillations and a very clear trend 
to the stabilization into a new stationary soliton, which 
does not belong to the family of the exact solutions. 



C. Collisions between embedded solitons 

Once knowing that the ELS of Eq. Q are stable 
against perturbations, it is natural to consider their ro- 
bustness against collisions. To this end, we simulated 
collisions between two ELS moving in opposite direc- 
tions. A typical example can be displayed for Eq = —4.8, 
p = 4.8, S = 10.3 and Ax = 0.4. In this case, all the 
solitons belonging to family 10-© have the same ampli- 
tude and width, Asoi = 0.38 and w — 1.5, but may differ 
in the velocity, which is controlled by the parameter Q. 
We can take a soliton moving to right, with Qi = 3.315 
and oi = 15.7, and a second one travelling to left, with 
Q2 = 7.9 and 02 = —9.0. Initially, these solitons are 
placed at nAx = —20 and nAx ~ 12, respectively. 

The collision between these solitons is shown in Fig. 
10, which demonstrates that they recover their shapes 
and velocities after the collision. In Fig. 11, we specially 
display the shape of the solitons after the collision (at 
t = 11). This figure shows that no radiation is emitted 
by the solitons after the collision, and it confirms that 
the solitons keep the same value of the amplitude, Agoi ~ 
0.38, which they had before the collision. At the moment 
shown {t = 11), the centers of the pulses are located at 
nAx — —78.6 and nAx = 139.2. Without the collision. 



they would be found (at the same moment of time) at 
the points nAx = —87 and nAx = 152.7, respectively. 
We thus conclude that the only effect of the collision is 
a shift in the position of the solitons (as occurs usually 
in integrable systems). The same result was observed in 
simulations performed at other values of the parameters. 



V. SUMMARY AND CONCLUSIONS 

In this work we have introduced a new one-dimensional 
dynamical lattice model. The underlying equation © 
includes next-nearest-neighbor couplings, and it may be 
regarded as a discretization of the complex modified KdV 
equation, or as an extension of the Ablowitz-Ladik model 
(in the staggered form). The model may apply to the de- 
scription of a Bose-Einstein condensate with dipole inter- 
actions between atoms, trapped in a deep optical lattice. 

A family of exact lattice-soliton solutions of the dis- 
crete equation was found, with the following characteris- 
tics: 

(i) These solitons can move with arbitrary constant veloc- 
ities (within a certain interval), without being hindered 
by the lattice discreteness. On the other hand, all the 
solitons have the same amplitude. 

(ii) All the soliton solutions are indeed embedded (with 
the exception of a part of the subfamily of zero-velocity 
solitons) because their internal frequency lies within the 
phonon band (calculated in a moving reference frame 
which moves along with the soliton, and contains the 
corresponding Doppler shift). 

(iii) The embedded solitons of Eq. Q are not isolated 
solutions (unlike most other models admitting embedded 
solitons), but form a continuous two-parameter family. 

(iv) Numerical simulations show that all the solitons are 
stable. 

(v) The exact soliton solutions with the fixed amplitude 
constitute a subset a larger family of stable discrete mov- 
ing solitons, which can be found in a numerical form. In 
particular, under a perturbation that reduces the soli- 
ton's norm, a perturbed soliton relaxes to a new one 
which belongs to the extended family. 

(vi) Collisions between the moving solitons known in the 
exact form are completely elastic (according to the results 
of systematic simulations) , leading solely to a shift of the 
solitons' positions. 

Some of these features, such as the existence of a con- 
tinuous family of moving lattice solitons and the elastic- 
ity of the collisions between them (as confirmed up to the 
numerical accuracy), suggest that the model may have a 
chance to be integrable. Investigation of this possibility 
poses a challenge for further studies of the model. 
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Figure captions 



The solitons moving to the right and left are determined 
by parameters, respectively, Qi — 3.315, oi — 15.7, and 
= 7.9, a2 = -9.0. 



Fig. 1. The phase velocity of the linear waves (times 
Ax) as a function of k [i.e., the right-hand side of Eq. 
(HH)] for £0 = 1 and Ax = 3. 



Fig. 11. The shape of the lattice field after the collision 



Fig. 2. The perturbed double-embedded lattice soli- shown in Fig 10 ai t = 11 
ton, at t = 20, which evolved from the initial condition 
((TB)l . with Aq = 0.64, w = 12.12, Q = 0.61, and Ax = 3. 
The parameters of Eq. (Q) are = I, p = —0.4, and 
5 — —0.85. In this figure and below, the shape of the 
soliton is shown as |r„| versus n. 

Fig. 3. The power spectrum of the numerical solution 
corresponding to the perturbed double-embedded lattice 
soliton shown in Fig. 2. 

Fig. 4. The phase velocity of the linear waves (times 
Ax) as a function of fc, for Eq = —4.8 and Ax = 0.4. 

Fig. 5. The perturbed single-embedded lattice soliton 
is shown aX t = 10. It has evolved from the initial con- 
dition lini, with Aq = 0.35, w = 1.503, Q = 7.9, and 
Ax = 0.4. The parameters of Eq. Q| are so = —4.8, 
p = 4.8, and S = 10.3. The velocity of the soliton is 
defined by the parameter a = 15.812. 

Fig. 6. The power spectrum of the numerical solution FIG. 1: 

corresponding to the perturbed single-embedded lattice 
soliton shown in Fig. 5. 

Fig. 7. Evolution of the amplitudes of two perturbed 
double-embedded lattice solitons. In both cases, the ini- 
tial conditions are given by Eq. H16|l . with w — 12.12, 
Q — 0.61, and Ax = 3. The parameters of Eq. are 
Eq — I, p — —0.4, and S — —0.85. The upper and lower 
curves correspond, respectively, to the initial amplitudes 
Aq = 0.64 and Aq = 0.48, which are larger and smaller 
than the value Agoi = 0.56 for the exact soliton solution 
in this case. 

Fig. 8. Evolution of amplitudes of two perturbed 
single-embedded lattice solitons. In both cases, the ini- 
tial conditions are taken as per Eq. H16|l . with w = 1.5, 
Q = 7.9, and Ax = 0.4. Parameters of Eq. (P) are 
Eq — —4.8, p = 4.8, and S = 10.3. The upper and lower 
curves correspond, respectively, to the initial amplitudes 
Aq = 0.415 and Aq = 0.345, which are larger and smaller 
that the value Agoi — 0.38 corresponding to the exact 
soliton in this case. 

Fig. 9. Evolution of the amplitude of a perturbed 
soliton corresponding to the parameters So — —I, p ~ 

0.4, S = 0.85 and Ax = 3. The initial pulse is take in the FIG. 2: 

form of Eq. ^ with Aq = 0.48, r = 0.61 and w = 12.12. 

Fig. 10. Collision of two single-embedded lattice soli- 
tons, at So = -4.8, p = 4.8, S = 10.3 and Ax = 0.4. 
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